﻿using System;
using System.Collections.Generic;
using DambachMulti;
using DambachMulti.DataStructures;
using LinearAlgebra;
using Microsoft.VisualStudio.TestTools.UnitTesting;

namespace MulticoreTests
{
    /// <summary>
    /// Summary description for LinearAlgebra
    /// </summary>
    [TestClass]
    public class UnitTests
    {        
        private TestContext testContextInstance;

        /// <summary>
        ///Gets or sets the test context which provides
        ///information about and functionality for the current test run.
        ///</summary>
        public TestContext TestContext
        {
            get { return testContextInstance; }
            set { testContextInstance = value; }
        }

        #region Additional test attributes
        //
        // You can use the following additional attributes as you write your tests:
        //
        // Use ClassInitialize to run code before running the first test in the class
        // [ClassInitialize()]
        // public   void MyClassInitialize(TestContext testContext) { }
        //
        // Use ClassCleanup to run code after all tests in a class have run
        // [ClassCleanup()]
        // public   void MyClassCleanup() { }
        //
        // Use TestInitialize to run code before running each test 
        // [TestInitialize()]
        // public void MyTestInitialize() { }
        //
        // Use TestCleanup to run code after each test has run
        // [TestCleanup()]
        // public void MyTestCleanup() { }
        //
        #endregion        

        [TestMethod]
		public void RunParallelGaussianInversionTest()
		{
			for( int i = 0; i < 100; i += 10 )
			{				
				//try {
				double[][] res1;
				double[][] res2;
				double[][] res3;
						
				Console.WriteLine("*** Gaussian Inversion ***: ");
							
				DateTime start;					
				double[][] Reg = ClassicMatrix.GetRandomMatrix( i, 36, 1, 1000 );
				start = DateTime.Now;	
				res1 = Invert( Reg );						
				Console.WriteLine( "Synchronous: ".PadRight(25) + DateTime.Now.Subtract( start ).TotalMilliseconds.ToString("f3").PadLeft( 8 ) + " ms" );					
				
				double[][] Par = ClassicMatrix.GetRandomMatrix( i, 36, 1, 1000 );
				start = DateTime.Now;
				res2 = InvertProcPooledWithStride( Par );						
				Console.WriteLine( "ProcPooledWithStride: ".PadRight(25) + DateTime.Now.Subtract( start ).TotalMilliseconds.ToString("f3").PadLeft( 8 ) + " ms");					
												
				double[][] Cur = ClassicMatrix.GetRandomMatrix( i, 36, 1, 1000 );
				Cur = ClassicMatrix.GetRandomMatrix( i, 36, 1, 1000 );
				start = DateTime.Now;
				res3 = InvertSmartPooled( Cur );						
				Console.WriteLine( "Dambach Multi-Smart: ".PadRight(25) + DateTime.Now.Subtract( start ).TotalMilliseconds.ToString("f3").PadLeft( 8 ) + " ms");					
					
				double diff = 0;
				if( ClassicMatrix.AreMatricesEqual( res1, res2, .1, out diff ) &&
					ClassicMatrix.AreMatricesEqual( res2, res3, .1, out diff ) ) //not efficient, but eh, who cares right?			
					Console.WriteLine( "Matrices(" + i + ") are equal! :)" );
				else					
					Console.WriteLine( "Matrices were not equal! :(" );
				Console.WriteLine( "Diff:" + diff );
				Console.WriteLine();				
			}
		}
		
		#region Inversions
		private double[][] Invert( double[][] A )
		{			
			double MatrixEpsilon = 0.001;	// value that we will count as zero
			int n = A.Length;				// the size of the square A
			double[][] B = ClassicMatrix.GetIdentityMatrix( n );

			for( int k = 0; k < n; k++) //for every row, starting with the first row(using 1 based indices!!!!!!!!!)
			{	//This loop does forward elimination			
				if ( Math.Abs( A[k][k]) < MatrixEpsilon ) // Check for small element on the diagonal, C# has A built in exception type 
					throw new DivideByZeroException( "Value less than epsilon(" + MatrixEpsilon +  ") on the diagonal, possible divide by zero condition! Cannot continue." );

				double xmult; 
				for( int i = k + 1; i < n; i++ ) {			//Get zeros for column k, //recurse over every row under the kth row
					xmult = -A[i][k] / A[k][k];				//what is the mult?
					ClassicMatrix.AddRowMultiple( A, B, xmult, k, i );	//reduce mat A, inverse mat B
				}
			}		

			double piv; //use just one reference to avoid excessive reference creation, less references = less garbage collection = faster code			
			for( int i = 0; i < n; i++ ) { //this loop performs normalization of the matrix...
				piv = A[i][i];
				for( int j = 0; j < n; j++ ) {
					A[i][j] /= piv;
					B[i][j] /= piv;
				}
			}

			for( int j = n-1; j >= 0; j-- ) //go backwards and solve for other x values			
				for( int i = j-1; i >= 0; i-- ) //I don't think we can parallelize this part as is				
					ClassicMatrix.AddColumnMultiple( A, B, -A[i][j], j, i );		
			return B;
		}
		
		private double[][] InvertProcPooledWithStride( double[][] A )
		{
			double MatrixEpsilon = 0.001;	// value that we will count as zero
			int n = A.Length;				// the size of the square A
			double[][] B = ClassicMatrix.GetIdentityMatrix( n );
		
			for( int k = 0; k < n; k++) //for every row, starting with the first row(using 1 based indices!!!!!!!!!)
			{	//This loop does forward elimination			
				if ( Math.Abs( A[k][k]) < MatrixEpsilon ) // Check for small element on the diagonal, C# has A built in exception type 
					throw new DivideByZeroException( "Value less than epsilon(" + MatrixEpsilon +  ") on the diagonal, possible divide by zero condition! Cannot continue." );								
				
				Multi.ProcPooledForWithStride( k+1, n, delegate( int i ) {			// Get zeros for column k
					double xmult = -A[i][k] / A[k][k];		// multiplier used on each row
					ClassicMatrix.AddRowMultiple( A, B, xmult, k, i );	// give our parallel code more to do
				}, 4 );
			}
	
			double piv; //use just one reference to avoid excessive reference creation
						//less references = less garbage collection = faster code			
			for( int i = 0; i < n; i++ ) { //this loop performs normalization of the matrix...
			//Multi.For( 0, n, 1, delegate( int i ) { //this loop performs normalization of the matrix...
				piv = A[i][i];
				for( int j = 0; j < n; j++ ) {
				//Multi.For( 0, n, 1, delegate( int j ) {
					A[i][j] /= piv;
					B[i][j] /= piv;
				} //);
			}
						
			for( int j = n-1; j >= 0; j-- ) //go backwards and solve for other x values
				for( int i = j-1; i >= 0; i-- ) //zeroing out whole columns of A at a time, and normalizing rows of B at the same time					
					ClassicMatrix.AddColumnMultiple( A, B, -A[i][j], j, i );
			return B;
		}
		private double[][] InvertSmartPooled( double[][] A )
		{
			double MatrixEpsilon = 0.001;	// value that we will count as zero
			int n = A.Length;				// the size of the square A
			double[][] B = ClassicMatrix.GetIdentityMatrix( n );
		
			for( int k = 0; k < n; k++) //for every row, starting with the first row(using 1 based indices!!!!!!!!!)
			{	//This loop does forward elimination			
				if ( Math.Abs( A[k][k]) < MatrixEpsilon ) // Check for small element on the diagonal, C# has A built in exception type 
					throw new DivideByZeroException( "Value less than epsilon(" + MatrixEpsilon +  ") on the diagonal, possible divide by zero condition! Cannot continue." );								
				
				Multi.For( k+1, n, delegate( int i ) {			// Get zeros for column k
					double xmult = -A[i][k] / A[k][k];		// multiplier used on each row
					ClassicMatrix.AddRowMultiple( A, B, xmult, k, i );	// give our parallel code more to do
				} );
			}
	
			double piv; //use just one reference to avoid excessive reference creation
						//less references = less garbage collection = faster code			
			for( int i = 0; i < n; i++ ) { //this loop performs normalization of the matrix...
			//Multi.For( 0, n, 1, delegate( int i ) { //this loop performs normalization of the matrix...
				piv = A[i][i];
				for( int j = 0; j < n; j++ ) {
				//Multi.For( 0, n, 1, delegate( int j ) {
					A[i][j] /= piv;
					B[i][j] /= piv;
				} //);
			}
						
			for( int j = n-1; j >= 0; j-- ) //go backwards and solve for other x values
				for( int i = j-1; i >= 0; i-- ) //zeroing out whole columns of A at a time, and normalizing rows of B at the same time					
					ClassicMatrix.AddColumnMultiple( A, B, -A[i][j], j, i );
			return B;
		}
		#endregion
	
		#region Matrix Basics 
        [TestMethod]
		public  void RunBasicMatricesTest()
		{
			double[][] A = ClassicMatrix.GetIdentityMatrix( 5 );			

			Multi.For( 0, A.Length, 1, delegate( int i ) {
				for( int j = 0; j < A[i].Length; j++ )									
        			A[i][j] += 1;
        	} );
			Console.WriteLine( "Added 1 to every cell in the identity Matrix stored result in A." );

            Multi.For(A.Length - 1, -1, -1, delegate(int i) {
                for (int j = 0; j < A[i].Length; j++)
                    A[i][j] -= 1;
            });
			Console.WriteLine( "Subtracted 1 to every cell in A, test was done in reverse Length-1 to 0." );
			ClassicMatrix.Write( A ); //visually compare the matrix to make sure it is the ident matrix
			double[][] B = ClassicMatrix.GetIdentityMatrix( 5 );
			double diff = 0;
			if( ClassicMatrix.AreMatricesEqual( A, B, double.Epsilon, out diff ) )			
				Console.WriteLine( "Matrices are equal. Test Passed." );
			else
				Console.WriteLine( "Matrices are NOT equal. Test Failed.");			

			for( int n = 1; n < 100; n++ )
			{
				double[][] one = ClassicMatrix.GetMatrix( n, (int i, int j) => Math.Sin(i) + Math.Cos( j ) ); //lambda expressions are teh hotness
				double[][] two = ClassicMatrix.GetZeroMatrix( n );

				Multi.For( 0, n, 1, (int i) => //  typing is necessary so we can prove things about the arguments!!!
					Multi.For( 0, n, 1, (int j ) => //AWESOMUH!  Not an efficiency test, but a correctness test
						two[i][j] = Math.Sin(i) + Math.Cos( j ) ) ); //now we compare the matrices
				if( ClassicMatrix.AreMatricesEqual( one, two, double.Epsilon, out diff ) ) //the max diff is less than the epsilon of the double, this is VERY acceptable
					Console.WriteLine( "Matsize(" + n + "):Sin(i) + Cos(j): Matrices are equal. Test Passed." );
				else			
					Console.WriteLine( "Matsize(" + n + "):Sin(i) + Cos(j): Matrices are NOT equal. Test Failed.");
			}
		}
		#endregion

        [TestMethod]
		public void RunDataStructuresTest()
		{
			Console.WriteLine( "Running Data Structures Tests" );
			const int testSize = 10000000;
			DateTime start = DateTime.Now;
			MultiArray<int> numList = new MultiArray<int>(testSize);
			for( var i = 0; i < testSize; i++ ) {
				numList.Add( i );
			}
			numList.Contains( testSize-1 );
			Console.WriteLine( "Multi: " + DateTime.Now.Subtract( start ).TotalMilliseconds );

			start = DateTime.Now;
			List<int> numList2 = new List<int>(testSize); //why is the implicit typing being touted as a feature!!!
			for( var i = 0; i < testSize; i++ ) {
				numList2.Add( i );
			}
			numList2.Contains( testSize-1 );
			Console.WriteLine( "Sync: " + DateTime.Now.Subtract( start ).TotalMilliseconds );			
		}
	}
}